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O Abstract 

^^ In this paper, we study the locaUzation phenomena in a slender cyhnder composed 

^ of an incompressible hyperelastic material subjected to axial tension. We aim to 

(-H construct the analytical solutions based on a three-dimensional setting and use 

Q^ the analytical results to describe the key features observed in the experiments by 

others. Using a novel approach of coupled series- asymptotic expansions, we derive 

^-H the normal form equation of the original governing nonlinear partial differential 

^ equations. By writing the normal form equation into a first-order dynamical system 

^vq and with the help of the phase plane, we manage to solve two boundary-value 

0^ problems analytically. The explicit solution expressions (in terms of integrals) are 

^ . obtained. By analyzing the solutions, we find that the width of the localization 

^""^ zone depends on the material parameters but remains almost unchanged for the 

r — same material in the post-peak region. Also, it is found that when the radius-length 

<^ ratio is relatively small there is a snap-back phenomenon. These results are well 

^ in agreement with the experimental observations. Through an energy analysis, we 

'kj^ also deduce the preferred configuration and give a prediction when a snap-through 

;_] can happen. Finally, based on the maximum-energy-distortion theory, an analytical 



criterion for the onset of material failure is provided. 
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1 Introduction 



The field of fracture meclianics is becoming extremely broad with the oc- 
currence of unexpected failure of weapons, buildings, bridges, ships, trains, 
airplanes, and various machines. There are two fundamental fracture criteri- 
ons: the strain energy release rate (i.e., G Theory) and the stress intensity 
factor (i.e., K Theory); see Arthur and Richard (2002). The experimentally 
determined stress intensity depends on the specimen size, and the fracture is 
accompanied by energy localization and concentration. 

Localization is manifested by degradation of material properties with localized 
large deformations, and this feature often results in formation and propagation 
of macrocracks through engineering structures. Due to the importance of lo- 
calization phenomena in structural safety assessment, much research has been 
conducted to resolve experimental, theoretical and computational issues asso- 
ciated with localization problems, as reviewed by Chen and Schreyer (1994) 
and Chen and Sulsky (1995). For hyperelastic materials, important progress 
has been made based on the gradient approach; see Aifantis (1984) and Tri- 
antafyllidis and Aifantis (1986). However, there is still a lack of analytical 
results for three-dimensional boundary-value problems. In this paper, hence, 
we study the localization in a slender cylinder composed of an incompressible 
hyperelastic material subjected to tension, based on an analytical approach 
to solve the three-dimensional governing equations. We also intend to provide 
mathematical descriptions for some interesting phenomena as observed in ex- 
periments. 

Jansen and Shah (1997) conducted careful experiments on concrete cylin- 
ders by using the feedback-control method. From two test series, the typical 
stress-displacement behavior for different height-diameter ratios with normal 
strength and high strength was obtained. It appears that the pre-peak seg- 
ment of the stress-displacement curves agrees well with the pre-peak part of the 
stress-strain curves, but the post-peak segment shows a strong dependence on 
the geometric size (i.e., the radius-length ratio). More specifically, the longer 
the specimen is, the steeper the post-peak part of the stress-displacement 
curves becomes. Also, they found that the width of the localization zone 
changes with the specimen size. In the experiment by Gopalaratnam and 
Shah (1985), it was found that the tangent value in the ascending part of 
the stress-strain curves seemed to be independent of the specimen size but in 
the post-peak part there was a softening region and no unique stress-strain 
relation. Schreyer and Chen (1986) studied the softening phenomena analyt- 
ically based on a one-dimensional model. Their results indicate that if the 
size of the softening zone is small enough (in a relative sense), the behavior 
of displacement-prescribed loading is unstable, and the softening curves are 
steeper than those with a larger size of the softening region. Here, we shall 
provide the three-dimensional analytical solutions to capture all the localized 
features mentioned above. 



Another purpose of this research is to provide a method judging the onset of 
failure in a slender cylinder subjected to tension. Here, we use the maximum- 
distortion-energy theory (the Huber-Hencky-Von Mises theory; see Riley et 
al. 2007), which depicts there are two portions of the strain energy intensity. 
One is the portion producing volume change which is ineffective in causing 
failure by yielding, and the other is that producing the change of shape which 
is completely responsible for the material failure by yielding. 
By constructing the analytical solutions for localizations, it is possible to get 
the point-wise energy distribution. Then, an expression for the maximum value 
of the strain energy can be obtained. With the Huber-Hencky-Von Mises the- 
ory, we can then establish an analytical criterion for identifying the onset of 
failure. 

Mathematically, to deduce the analytical solutions for localizations in a three- 
dimensional setting is a very difficult task. One needs to deal with coupled 
nonlinear partial differential equations together with complicated boundary 
conditions. Further, the existence of multiple solutions (corresponding to no 
unique stress-strain relation) makes the problem even harder to solve. Here, 
the analysis is carried out by a novel method developed earlier (Dai and Huo, 
2002; Dai and Fan, 2004; Dai and Cai, 2006), which is capable of treating the 
bifurcations of nonlinear partial differential equations. Our results yield the 
analytical forms of the strain and stress fields, total elongation, the potential 
energy distribution and the strain energy distribution, which are characterized 
by localization phenomena. In particular, it is found that once the localiza- 
tion is formed its width does not change with the further increase of the total 
elongation, which is in agreement with the experimental observations. We also 
provide a description for the snap-through phenomenon. 
The remaining parts of the paper are arranged as follows. In section 2, we 
formulate the three-dimensional govering equations for the axisymmetric de- 
formations of a circular cylinder. We nondimensionalize them in section 3 to 
identify the key small variable and key small parameters. Then, in section 4, 
a novel method of coupled series-asymptotic expansions is used to derive the 
normal form equation of the original system. By the variational principle, in 
section 5, we derive the same equation by considering the energy. In section 
6, with the help of the phase plane, we solve the boundary-value problems for 
a given external axial force and a given elongation, respectively. In section 7, 
through an energy analysis, we determine the most preferred configurations 
and give a description of the snap-through phenomenon. Also, an analytical 
criterion for identifying material failure based on the Huber-Hencky-Von Mises 
theory is discussed. Finally, concluding remarks and future tasks are given in 
section 8. 



2 Three-dimensional governing equations 



We consider the axisymmetric deformations of a slender cylinder subjected to 
a static axial force at one plane end with the other plane end fixed. The radius 
of the cylinder is a and the length is /. We take a cylindrical coordinate system, 
and denote (i?, G, Z) and (r, 9, z) as the coordinates of a material point of the 
cylinder in the reference and current configurations, respectively. The radial 
and axial displacements can be written as 



u{R, Z) = r{R, Z) - R, w{R, Z) = z{R, Z) - Z. 



(1) 



We suppose that the cylinder is composed of an incompressible hyperelastic 
material, for which the strain energy density $ is a function of the first two 
invariants Ji and I2 of the left Cauchy-Green strain tensor, i.e., $ = <l>(/i,/2). 
Moreover the first Piola-Kirchhoff stress tensor S is given by 



(9$ 



-pF- 



(2) 



where F is the deformation gradient and p is the indeterminate pressure. If 
the strains are small, it is possible to expand the first Piola-Kirchhoff stress 
components in terms of the strains up to any order. The expressions for the 
stress components are very lengthy, and due to the complexity of calculations, 
we shall only work up to the third-order material nonlinearity. The formula 
containing terms up to the third-order material nonlinearity has been provided 
by Fu and Ogden (1999): 
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Zjjj — Cijii}J]kl + 2^ji!-knmVklVmn + r,(^jilknmqpVklVmn'npq 



(3) 



+Po{Vji - VjkVki + VjkVklVli) - P*{^ji - Vji + VjkVki) + 0{\r]ij\'^), 



where 
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Po is the pressure value in Eq. (2) corresponding to zero strains, p* is the 



incremental pressure, and a] 



a 



'jilky "^jilknm^ 



and a?, 



jilknmqp 



are incremental elastic 



moduli defined by 



1 _ _0^f_| 2 _ 8^$ I 



,3 _ y$ 

^ijdFj^ldFmnOPpq ' 



(^jilknmqp dFijdFkidF,„ndFpr l^=-f- 



It can be found tliat 

po = 4$oi + 2$io, 
wliere $oi denotes tlie first-order partial derivative of $ with respect to the 
invariant /2 at F = I, $io denotes the first-order partial derivative of $ with 
respect to the invariant Ji at F = I. In the following derivations, we shall also 
use $jj to denote the i-th order and the j-th order partial derivative of $ with 
respect to the invariants /i and /2 at F = I. All the coefficients in Eq. Q can 
be expressed in terms of $io, $20, "^oi, "^02, and $11, and here for brevity the 
expressions are omitted. 
The equilibrium equations for a static and axisymmetric problem are given by 

The incompressibihty condition yields that 

Vii = ^VmnVnm " ^vl " det(r7„„). (7) 

We consider the case where the lateral surface of the cylinder is traction-free, 
then the stress components T^rR and S^z should vanish on the lateral surface. 
Thus, we have the boundary conditions: 



^rR\R=a 



0, Kz\R=a = 0. (8) 



Eqs. (l5])-([7]) together with Eq. ^ provide three governing equations for three 
unknowns u,w and p*. The former two are very complicated nonlinear par- 
tial differential equations (PDE's) and the boundary conditions (tsl) are also 



complicated nonlinear relations (cf. ([ll|)-(15)). To describe the localization, 
one needs to study the bifurcation of this complicated system of nonlinear 
PDE's. As far as we know, there is no available mathematical method. Here, 
we shall adapt a novel approach involving coupled series-asymptotic expan- 
sions to tackle this bifurcation problem. A similar methodology has been de- 
veloped to study nonlinear waves and phase transitions (Dai and Huo, 2002; 
Dai and Fan, 2004; Dai and Cai, 2006). First, we shall nondimensionalize this 
system to identify the relevant small variable and small parameters. 



3 Non-dimensionalized Equations 



We introduce the dimensionless quantities through the following scales: 



h V* 

s = I'^s, Z = Iz, w = hw, V = —V, — 



I 



P , 



(9) 



where / is the length of the cylinder, h is a characteristic axial displacement 
and /i is the material shear modulus, with a transformation being defined by 



u = vR, s = R . 



(10) 



Substituting Eqs. (g, Q and (|lO|) into Eqs. <Q-^, we obtain 

-pl + buv^ + ibioWs + h2Wzz + s{hiAVsz + '^hi^Wss) H = 0, 

- 2p* + %h2Vs + hi^Vzz + hiiWsz + sAh2Vss H = 0, 

2v + Wz + 2sVs + e[v'^ + 2vWz + 2s{vvs — VzWg 

+VsWz)] + e^[v'^Wz + 2sv{VsWz - VzWs)] = 0, 



(11) 

(12) 
(13) 



h 
where e = — is regarded as a small parameter. For convenience, we have re- 
placed s, V, w, z, p* by s, V, w, z, p* in the non-dimensionalized equations. Here 
and thereafter, a subscript letter is used to represent the corresponding partial 



derivative (i.e., Vz = |^). The full forms of (11) and (12) are very lengthy and 
we do not write out the nonlinear terms explicitly for brevity. 
Substituting ^ and (10) into the traction-free boundary conditions ([s]), we 
have 

—p* + 6i4f + 2biWz + i'2h2Vs + e[p*v + &46f ^ + hgVWz + h^vol 

+u{2p*Vs + hi^vVs + 364f^ — 2buVzWs + I2biw'^ + Ab^VsiUz) 

+iy'^4b5vl] + e^[-p*v'^ + 3b47V^ + 30bjv^Wz + bi^vwl + 463 wf 

+ l'{-2p*{2vVs + VzWs) + b22v'^Vs + Ibgwl + b23VVzWs (14) 

+2SbQVw'^ + 72bi2VVsWz + 2bigv1Wz + 2b^VzWaWz + Sfoigwf w 



s"^z 



+bi^Vswl) + u'^{-Ap*vl + ?>b^ivvl + Abi^Vsvl + 9,bi^VsVzWs 
+l<obisiVswl + b2ovlwz) + v-^b2ivl]\s=u = 0, 
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+12b4WsWz + iy{-2biiVsV;, + 2Ab4VsWs)] + e'^[-p*{vv;, + v^Wz) 

+bMv'^Vz + 2biQV^Ws + b^QVV;,Wz + b32VWsWz + &l3f ^w^ 
+b33WsW'l + v{-2p*VsVz + b23VVsVz - b^vl + ^Qb^VVsWs 

+2biivlws - l2b3Vzwl + IQb^wl + 2bsVzVsWz 



+16biQVsWsWz) + ij'^{Abi3vlvz + 166i4ffw 



^SJ\ \S = U ") 



0, 



(15) 



where z/ = ya- is a small parameter for a slender cylinder. 



Then, Eqs. (11)-(15) compose a new system of complicated nonlinear PDE's 



with complicated boundary conditions, which is still very difficult to solve 
exactly. However, it is characterized by a small variable s and two small pa- 
rameters {e and z/), which permit us to use expansion methods to proceed 
further. 



Remark: The coefficients bi,b2,--- in Eqs. (11)-(15) can be expressed in 



terms of $10, $20, "^oij "^02 and $11, and for brevity we omit their expressions. 



4 Coupled Series- Asymptotic Expansions 



We note that s is also a small variable as < s < u. An important feature of 



the system (11)-(15) is that the unknowns w,v, and p* become the functions 



of the variable z, the small variable s and the small parameters e and u, i.e.. 



w = w(z, s; e, z/), f = t>(z, s; e, z/), p* = p*{z,s;e,h'). 



(16) 



To go further, we assume that w, v,p* have the following Taylor expansions in 
the neighborhood of the small variable s = 0: 



w 



p* = Po{z; e, u) + sPi{z; e, v) + s^P2{z; e, u) + 
V = Vo{z; e, z/) + sVi{z; e, z/) + s'^V2{z; e, u) + 

2t 



Wo{z; e, v) + sWi{z] e, u) + 3^2(2;; 5, z/) + 



(17) 
(18) 
(19) 



Substituting Eqs. (17)-(19) into Eq. (12) and equating the coefficient of s° to 
be zero yields that 

-2Pi + 862 Vi + bioVo,, + hi^Wi, + £(2PiK) + 2P^zWi + +2Po(4\/i + W^i.) 

+hsWl - ^h^iVoMi + h^Vl + SfogV^iW^o. + 664K)..W^o. + 26nW^o..W^i + 664^.1^0.. (20) 



Similarly, substituting Eqs. (17)-(19) into Eq. (11) and setting the coefficients 



of s° and s^ to be zero, we obtain 



^-2h^VoWo,, + 2464^^1 W^o. + &2oK).W^o. + 2627^^0.1^0..) + e^H^iVo, Vi, W^o, W^i, A) = 0, 



(21) 



-Pi, + 166iol^2 + 2614V1, + 62^^!.. + e{Pi,Wo, + Pi(4K). + Wo,,) + Po.W^i. 
+Po(4\/i, + W^i,,) - 26nW^iV^o.. + 664^).^).. + 96641^2^^0. + 262oVl.W^o. + 7264^^1 W^i.) 
+V^{lQhW^ - be^Vo, + SbeWo,,) - be^Vo.Wi, + 263(^^0.1^1.). 
+Vo{16buW2 + 2626^. + 4661^1,,) + e^HM, V,, Wo, W,, Pq, Pi) = 0. 



(22) 



The expressions of Hi , H2 and H^, are very lengthy, whose expressions are 



omitted for brevity. From the in compressibility condition (13), the vanishing 



of the coefficients of s" and s^ leads to the following two equations: 
2H + 1^0. + ^K)(K) + 2Wo,) + eX'Woz = 0, 
m + Wu + e{-2WiVo, + AVM + Wo,) + IVoW^,) 
^e-'Vo^-miVo, + WiWo, + VoWi,) = 0. 



(23) 



(24) 



Substituting Eqs. (17)-(19) into the traction-free boundary conditions (14) 
and (15), and neglecting the terms higher than 0(z/£,£^), we obtain 



-Po + buVo + 261W0. + i^i-Pi + 67014 + 26iW^] 



i.J 



-e{PoVo + 646^0' + &39l^ol^o. + 66ll^( 



Oz) 



(25) 



+e^H,{Vo, Wo, Po) + lyeH.iVo, Vi, Wo, Wi, Po, Pi) = 0, 



ho{2Wi + Vo,) + ubioiVu + 41^2) + eiPoVo, + Vo(-64iK). 

+2h7Wi) + Wo.i-bnVo. + 1264^^1)) + e^He{Vo, Wo, Wi, Pq) (26) 

+iyeHr{Vo, Vi, Wo, W,, W2, Po, Pi) = 0, 



where the lengthy expressions for 7^4 — if 7 are omitted for brevity. Eqs. (20)- 



(26) are seven nonhnear ordinary differential equations, which are the govern- 



ing equations for the seven unknowns Pq, Pi,Vo,Vi,Wo,Wi and W2- Mathe- 
matically, it is still very difficult to solve them directly. To go further, we shall 



use the smallness of the parameter e. From Eq. (23), we obtain 



Wf 



Oz 



-K)(2 + £Vo)/(l + £K)) 



-2K) + SVffe - W^e' 



(27) 



Using the above equation in Eq. (24), we can express Vi in terms of Vq and 



Wi, and then from Eqs. (20) and (22), we can express Pi and W2 in terms 



of Vo, Wi and Pq. Substituting the expressions for Wo, V\, P\ and W2 into Eq. 



(21), Poz can be expressed in terms of Vq and W\. Then, from Eq. (25), P( 



can be expressed in terms of Vq and W\. Finally, from Eqs. (21) and (|26|), we 
obtain 



46io(W^i - V^o.) + eVo-f{Woz - 2Wr) + z/^(K),,, + W^,,) + e'^V^iho^Voz - 692^^1) 



(2^ 



nos- 



+i/£(8627W^iW^u + h2^Voz{Vozz + 2Wu) - &ioH^iK).. + &29K)W^i.. + ^VoVozzz) = 0, 



+ve{-^W^W^z - hoVozW,, + -hoWiVozz - fesiK).^).. - &32K)W^i.. + ^K)K)...) = 0. 

o Z o 



(29) 



We note that the above two equations come from the axial equilibrium equa- 
tion (the coefficient of s") and the zero tangential force at the lateral surface, 
the two most important relations for tension problems in a slender cylinder. 



By eliminating Wi^z from Eqs. (28) and (29) and then expressing Wi in terms 



of Vo, finally we obtain an equation for the single unknown Vq: 



3&10M). - 3eb28VoVoz + Se^feiioT/gSy^^ _ ^jjb^^Vozzz + yeh^iyozVozz + VoVozzz) = 0. (30) 



By further using Eq. (27), we obtain the following equation for the axial strain 
Wo,: 



36ioW^o.. - ehio^WoMozz + 'ie%uoWlW^ 



Oz^^Ozz 



-ji^hoWozzzz + i^eb35{4Wo,,Wozzz + ^WozWq,,,,) = 0. 



(31) 



Integrating Eq. (31) once, we obtain 



1 Q 

ShoWoz - e^Wl + e%noWi^, - -^vhoWozzz + ^eh.^Wl, + 2WozWozzz) = C, (32) 



where C is an integration constant. To find the physical meaning of C, we 
consider the resultant force T acting on the material cross section that is pla- 
nar and perpendicular to the cylinder axis in the reference configuration, and 
the formula is 



2-K a 



T= f f J^zzRdRde. 



(33) 







After expressing E^^ in terms of Woz by using the results obtained above, the 
integration can be carried out, and as a result we find that 



T = Sna^i^eiSbioWo, - e'-f^W^ + e%ioW^ 

-l^W^^Vozzz + ^^£^35(^^0.. + '^WozWozzz)). 



3 
Oz 



(34) 



Eq. (34) as 



Comparing Eqs. (32) and (34), we have C 



T 



Sna'^He 



Thus, we can rewrite 



3610^1^0. - "TieWoz)^ + buoieWozf - ^i^eb^oW, 



Ozzz 



(35) 



+ ^.£^635(^^0'.. + 2WozWozzz) = ^. 



If we retain the original dimensional variable and let V = Wqz = bWqz, we 
have 



V + D,V^ + D2V^ + a\-^-Vzz + D,yl + 2D^VVzz) = 7, 



(36) 
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where 



Di 



noo 



66 



D, 



10 



h 
36 



10 



D, 



^35 



10 



36 



-,7 



10 



T 

246iovra2/i 



(37) 



Since Eq. (36) is derived from the three-dimensional field equations, once 



its solution is found, the three-dimensional strain and stress fields can also 
be found. Also, it contains all the required terms to yield the leading-order 



behavior of the original system. Therefore, we refer Eq. ( 36 ) as the normal form 



equation of the system of nonlinear PDE's (11)-(13) together with boundary 



conditions (14) and (15) under a given axial resultant. 



5 The Euler-Lagrange Equation 



It is also possible to deduce the equation for V = sWqz by considering the 
total potential energy and then using the variational principle. By using the 
expansions obtained in section 4, we can express the two principal invariants 
/i and I2 in terms of Wqz- The results are 



j^ _ 3 = ^,^wl - 2e'Wl + 2e'W^^ + s[e\^^^Wlz - ^^W^-M.^..) 



+e 



3(13$oi + 22<l> 



10, 



»^v^^.+ ^'!;M"':°W '^w)i. 



32($oi + $io) ""'^ 8(<l>oi + $io 

J2 - 3 = ?>e^Wl - Ae^'Wl, + ^e^W^, + s[e 



2w2 ..3w3 , ,.4w4 , .^.2rlw2_^A^^^^^ 

64 F 



ozyyozzz) 



^sH- "^'^^^'l'"^ W,.WL + '%*51±|^„.. „,„.)]. 



64($oi + $10) 



8($oi + <^ 



lOj 



(3J 



From Eqs. (38 ), we know the first terms in the right-hand sides of Eqs. (38 ) are 



second-order nonlinear. To be consistent with the third-order material nonlin- 
earity of the stress components, the strain energy $ should be expanded up 
to the fourth-order nonlinear terms. So, according to the Taylor's expansion, 
we have 



$ = $io(/l - 3) + $01 (/2 - 3) + i[(/i - 3)2$20 

+2(/i - 3)(/2 - 3)<l>ii + ih - 3)2$o2] + ■ ■ ■ . 



(39) 



11 



In Eq. (39) it is sufficient to keep the second-order terms of /i — 3 and I2 — 3. 



Substituting the expressions of Ji and I2 in Eqs. (|38|) into Eq. (39), we have 

(40) 



2 
3 



$ = e^^i{Qh,^Wl - ^-ehr^oWl + e'^no^J + se^ iiQ^W.W^ 



'64 



15 
64 



&ioW^o.W^o... + e{F^WoMl^ + FsW^o'.W^o...)]. 



The stored energy per unit length is given by 

a 27r 

^ = f f ^RdRdQ. 





(41) 



Substituting Eq. (40) into Eq. (41) and carrying out the integration, we obtain 



the average stored energy over a cross section: 



81. 



^ = ^ = 2e^fi{6hoWl - -ebiooW^, + e%noW^,) + 2aV/i[-6ioW^( 



64 



Ozz 



-]^hoWo,Wo,., + £(FiWo.W^o.. + F2WIW0,,,)]. 
64 



(42) 



Letting V = Wqz = ^Wqz, we can further rewrite the above equation as 



27 
vl> = E[lV' + \D^V^ + \D2V^ + a\^Vl 



256 ^ 256 



VVzz + F^VVl + F2V^Vzz)l (43) 



where E = 24/i6io is the Young's modulus, Fi and F2 are constants related to 

material parameters. 

The total potential energy for a force-controlled problem is given by 



L = Tra\ I ^dZ -E -fVdZ) 





na'EJi-^V + \V' + \d,V' + \d2V' + a\^^V, 



(44) 



-iaVVzz + F^VVi + F2l^Vzz))rf^. 



By the variational principle, we have the following Euler-Lagrange equation: 
dL d dL d^ dL 



+ 



dV dZ dVz dZ'^ dVzz 



0, 



(45) 



12 



which yields that 



V + D^V^ + D^V" + a\D^Vl - -Vzz + 2D^VVzz) = 7, 



(46) 



which is just Eq. (36). This shows that the normal form equation (36) obeys 



the variational principle for energy. 

If we multiply Vz to both sides of Eq. (46), it can be integrated once to yield 

that 



1^.2 1 



1 

4^ 



,2/-'-ta2 



^V + ^DiV^ + ^D2V^ - a\^Vl - D^VV^) = ^V + K, 



(47) 



where K is an integration constant. 

In the following section, we shall discuss the solutions for two boundary-value 



problems based on Eqs. (46) and (47), and reveal their main characteristics. 



6 Solutions for two boundary- value problems 



We rewrite Eq.(46) as a first-order system as follows: 



Vz = y, 



(4J 



yz 



a\\-2D,V) 



Without loss of generality, we take the length / of the cylinder to be 1, then a 
is equivalent to the radius-length ratio. We suppose that the two plane ends 
of the cylinder are attached to rigid bodies. Then we have 



z = {or constant), 



at Z = 0, 1, 



(49) 



and 



R, 



at Z = 0, 1. 



(50) 



We point out that although Eq. (46) is one-dimensional, it is derived from the 



three-dimensional governing equations, and as a result we can also derive the 
proper boundary conditions by considering the condition in the other (radial) 



dimension such as Eq. (50 ). If one directly introduces a one-dimensional model 



(say, using a gradient theory), such an option is not available. So, this is 
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another advantage of Eq. (46). 



From Eqs. (49) and (50), we have 



wr = (i.e., zr = 0), and ur = (i.e., vr = 1) at Z = 0, 1. 



(51) 



Substituting Eq. (51) into Eq. (JTJ) and integrating with respect to R once, we 
obtain 



Wz =0, 



at Z = 0,1. 



(52) 



Thus, the proper boundary conditions for Eq. (46) are 
V" = 0, at Z = 0A. 



(53) 



To solve this boundary- value problem of the first-order system (48) under Eq, 



(53), we shall conduct a phase-plane analysis with the engineering stress as 



the bifurcation parameter. The critical points of this system is given by 



y = 0, and V + DiV^ + D2V^ 



1- 



(54) 



Here we shall consider a class of strain energy functions $(Ji,/2) such that 



the 7 — V^ plot based on Eq. (54)2 has one peak, and this requires that 



Di < 0, D2> 0, Dl > 4D2. 



The J — V curves corresponding to Eq. (54)2 are plotted in Fig. 1. 




Fig.l. y-Vplot 



In this figure, 7^ is the local maximum of the stress and Vm is the corre- 
sponding strain value, and they are given by Kn = 3^,^ -, 7^ = 

_ 2D^+2(D^-3D^f^^-9DW, ^ ^^^^ ^^ ^^^^ ^^ ^ -9.45,1^2 = 22 and D, = -2, 

Vm = 0.07; when we take Di = —6.65, D2 = 11 and D3 = —2, Vm = 0.1; when 
we take Di = —5.53, D2 = 7.6 and D^, = —2, Vm = 0.12. The three curves in 
Fig. 1 correspond to these values of Di, D2 and D3, respectively. 
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In the following discussions we consider the tension case so that 7 > 0. Similar 
analysis can be made for the compression case, which will not be discussed 



here. Equation (47) can be rewritten as 



^ Z 



a%-D,V) 



(55) 



In this paper, we consider the case of D3 < 0. New phenomena can arise for 
1^3 > and the results will be reported elsewhere. For the present case, the 
phase plane always has a saddle point and a center point as 7 varies, which is 
shown in Fig. 2. 




Fig. 2. The phase plane 



In this figure, (K, 0) and (Vc, 0) are a saddle point and a center point, respec- 
tively. There are two solutions for the same stress 7, which are represented 
by the curve 1 and the curve 2 in Fig. 2, respectively. For curve 1, the right 



hand of Eq. (55) have four real roots, which we label in an increasing order 



by ai, (71, (72 and 0^2. We note that the smallest root cti is smaller than Vs. So, 



from Eq. (55) we obtain the following expression: 



Vz 



V2D'2 



Wi - 8D3V 



^{V - a,){V - g,){V - g2){V - a^). 



(56) 



Then, an integration leads to 



Z = < 



ai 



\., , 57 

1 a r / 1 - SDat 1 

2 + v^^i V(i-«i)(i-^i)(i-^2)(t-«2) ' ^^2'^J- 
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By Eq. (53), ai can be determined by the following two equations: 



1 
2 



ai 



-D. 



8D3 



\\-K-^t+ \t^ + \D^t^ + \D2t* '' 



(5J 



K = -701 + -al + -Dial + -D2a\. 



(59) 



Once ai is found, the solution corresponding to curve 1 can be obtained from 



Eq. (57) by numerical integration. In Fig. 3, we have plotted the solution 



curves for three different values of the engineering stress 7. 




Fig. 3. V-Z plot. From bottom to top, 
y=0. 037825, y=0. 04005, y=0.0445 



From this figure, we find there is nearly a uniform extension in the middle part, 
but there are two boundary-layer regions near the two ends of the cylinder in 
order to satisfy the boundary conditions. 
There is another solution which is represented by curve 2 in Fig. 2, and we 



denote the point as a at which Vz = 0. Then, Eq. (55) can be rewritten as 



Vz 



V^D^ 



'aVl-SD^V 



^ia~V)i(3-V)[iV-m)^ + n% 



(60) 



where f3 is another real root of the right-hand of Eq. ( 55 ) , and m 



4:Di+3{a+l3)D2 



2 _ -WDj + 24(a + P)DiD2 + 9^2(8 + (Sa^ + 2ap + 3p^)D2) 



6D2 



and n 

Then, we obtain 



36D^ 



Z = I 



1 ^ f l-SDst 1 

2" v^i^/ V (a - t)(/3 - t) [{t - m)2 + n2] ' ^ ^ ' 2^ 



1 a 

7^ + 



2 ' ^/2D~2J \ {a-t){P-t)[{t-my + n^] ^' ^ ^2'^^' 



(61) 
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By Eq. (53), a can be determined by 



1 
2 " 

K 



-D, 



SD'i 



1„ , 1 



-7a + -a^ + -Dia^ + -D2a^. 



(62) 
(63) 



By numerical integration, we can get a from Eqs. (63) and (62). Then the 



solution corresponding to curve 2 can be obtained from (61) by numerical 



integration. In Fig. 4, we have plotted the solution curves for three different 
values of the engineering stress 7. 




Fig. 4. V-Z plot. From top to bottom (along the dashed line) 
7=0.037825, y=0. 04005, y=0.0445 

In this figure, there is a sharp-change region in the middle of the slender cylin- 
der, that represents the localization and energy concentration. Moreover, the 
tip is sharper when the engineering stress is smaller. 



From Eq. (61 ), one can see that the localization solution depends on Z through 
the form (Z — ^)/a, and this implies that the localization zone width is pro- 
portional to a for a fixed length; see Jansen and Shah (1997). 
The solutions obtained above are for a given 7. To obtain the solutions for a 
displacement-controlled problem, we follow the idea in Dai and Bi (2006). For 
that purpose, we need to get the engineering stress-strain curve. 
The total elongation is given by 



Wo 



z=i 



Wo 



z=o 



VdZ = A, 



(64) 



where the total elongation A is actually the engineering strain since we have 
taken the length of the cylinder to be 1. According to the symmetric phase 



plane and Eqs. (57) and (61), V^ is a function of Z, so we can get the total 



elongation by numerical integrations. In Fig. 5, we have plotted the curves 
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between the total elongation A and the engineering stress 7 with different 
material coefficients corresponding to Fig. 1. 




- V,„=0.07 

- K„=0.1 
V.=0.12 



Fig. 5. The engineering stress-strain curves 



Segment 1 corresponds to the solution given by Eq. (57) (we call it as Solution 
1), and segment 2 corresponds to the solution given by Eq. ( [6l] ) (we call it as 
Solution 2). For a displacement-controlled problem (i.e., given A), from Fig. 
5, we can get the corresponding 7 value(s), then the solution(s) is given by 



Eq. (|57|) or Eq. J6l|. 

From Eqs. (1), ([9j), (10), (18) and (19), we can get the shapes of the cylinder 
corresponding to Eq. (61) under different material coefficients for a given A, 
which are shown in Fig. 6, where we take D^ = —0.5, F2 = —4, and a = 0.04. 




7=0.023146, 6=0.216 y=0. 0340674, 5=0.22 

(a) A=0. 0740121 for Vm=0.07 




7=0.0281586, 5=0.25 

(b) A=0. 081 1761 for Vm=0.1 



7=0.0356, 5=0.236 



7=0.0325236, 5=0.24 7=0.0454186, 5=0.24 

(c) A=0. 10055 forVm=0.12 

Fig. 6. Shapes of the cylinder corresponding to the localization 
solutions and the comparison of the localization widths 5 



In this figure, the width 6 of the locahzation is defined as 1 — 2Z, where Z 
is the point where the rate of the slope of the surface radial displacement is 
the maximum. From the above figure, one can see that for different material 
coefficients the localization widths are different and the localization width 
is almost the same for the same material coefficients with different loads of 
engineering stress. That is to say, for different materials the localizations have 
different widths, but for the same material, the localization width is almost 
uniform during the loading process. 

Here and thereafter, we fix the material constants to be Di = —6.65, D2 = 11, 
and D^ = —2. By the same way, we can get the relations between the total 
elongation A and the engineering stress 7 with different radius-length ratios, 
which are shown in Fig. 7. 
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Fig. 7. The engineering stress-strain curves for different a 



In this figure, we observe that there is a snap-back for a relatively small value 
of a. We also see that the point c (across which there are multiple values 
for 7 for a given A) moves up and toward right as the value of a increases. 
For example, when a = 0.03, •jc = 0.0339, and when a = 0.04, 7c = 0.0399. 
The post-peak curves show very significant changes. There is no unique stress- 
displacement relationship in the post-peak region. The thinner the specimen 
is, the steeper the curve becomes, which is in agreement with the experimental 
results by Jansen and Shah (1997). From this figure, we see that for a = 0.03, 
unstable behavior is predicted for a displacement-controlled loading whereas 
larger values of a yield results that are stable. Similar conclusions follow from 
the examples given by Schreyer and Chen (1986). 

As to a = 0.045, there is a stable relation between the total elongation A and 
the engineering stress 7, which is in agreement with the experiment result by 
Gopalaratnam and Shah (1985), who conducted tensile tests on concrete under 
carefully controlled loading conditions and with refined measuring techniques. 
We note that for a displacement-controlled problem, after the elongation A > 
Ac (cf. Fig. 7) there are bifurcations from one solution to two solutions (at 
A = Ac), to three solutions (Ac < A < Ap), to two solutions (A = Ap), and 
to one solution (A > Ap). The shapes of the cylinder corresponding to these 
solutions are shown in Fig. 8 for F2 = —4, and a = 0.03. 
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y=0.0445465 



y=0.0302207 



7=0.037825 



(b) Ac<A=0.086133<Ap 



7=0.0439561 



(c) A=Ap=0.0927046 



7=0.0259625 



7=0.0259623 



(d) A=0.0927053>Ap 
Fig. 8. The shapes of the cylinder under different elongations 



The above figure also manifests tliat tlie middle of the cylinder becomes thin- 
ner than the two ends as we pull the slender cylinder. The middle part is 
thinner as the engineering stress decreases for a given A, which agrees well 
with the experimental results. 
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7 Energy analysis and failure criterion 



As discussed in the previous section, for a relatively small a, there are multiple 
solutions for A > Ac. Of course, in reality only one solution can be observed 
at one instant. In this section, we shall further consider the energy values for 
these solutions to deduce which one is most preferred. 



From Eq. (|47|), we have 

7 
3a2 (-1 + 8^3'^) 



^. 2{l2K + l2-iV-W^-ADiV''-?,D2V^) 

^z = 0^9/ 1 , on TM • (65) 



Substituting Eq. (65) into Eq. (46), we obtain 



^^^ = - o/ i^tn 1/^2 (37 + 247^^3 - 3F - W^V^ 

3(-l + 8/^3^)^ (66) 



Then by using Eq. (44), we can express the potential energy (for a given 7) 
in terms of V (scaled by na^E): 



384(-l + 8^3^)2 

-256(4 + 3a^)F2V + 16L'3V'(123 - 5a'^ + 256(2 + a^)F2V)) 

+3{D2V\~59 + lOa^ - 4096^1^2 - 512(1 + a^)F2V 

+4DsV{2AQ - 15a2 + 256(4 + 30^)^2 \/)) + 2(54fs: - 59V^ ^57) 

+5aV2 + II87I/ - 5a^-fV + 4096V\K - V{V - 2-f))Dl 

+25QV{AK + V{~{2 + a^)V + (4 + a^h))F2 

+4D3\/(-2(118 + 5a^)K + V{2ASV - ha^V - 4927) 

+256l^(2(-4 + a')K + F((4 + a')V - 87))F2)))). 

The stored energy is given by 

G = P^-iV. (68) 



Then from Eqs. (58), (61), (67) and (68), one can calculate the energy dis- 
tributions for a given elongation. The stored energy curves corresponding to 
those values of A in Fig. 8 are plotted in Fig. 9. In this figure, labels 1, 2 and 
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3 correspond to different values of 7 (in a decreasing order). 




(a) A=Ac^0.0849051 




(b) Ac<A=0.086133<Ap 




(c) A=Ap=0.0927046 




(d) A=0.0927053>Ap 
Fig. 9. The stored energy distribution curves 



In Fig.9(a), the total stored energy values for curve 1 and curve 2 are re- 
spectively G\ = 0.00241013 and G^ = 0.00206265. Thus for a displacement- 
controlled problem, as G2 < G\, the shape in the right of Fig. 8(a) represents 
a preferred configuration, and at A = Ac there could be a bifurcation from 
Solution 1 to Solution 2 (a localization solution). Correspondingly, there is a 
snap-through in the engineering stress-strain curve at A = A^. 
In Fig.9(b), the total stored energy values for curve 1, curve 2 and curve 3 



are respectively G\ = 0.00247747, 



G2 



0.00230487 and G* 



0.00186959. 
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For a displacement-controlled problem, as G\ is the smallest, the shape in the 
bottom of Fig. 8(b) represents a preferred configuration. 
In Fig.9(c), the total stored energy value for curve 1 and curve 2 are respec- 
tively G\ = 0.00275441 and G^ = 0.00170114. For a displacement-controlled 
problem, as G2 < G\, the shape in the right of Fig. 8(c) represents a preferred 
configuration. 

In Fig. 10, we have plotted the engineering stress-strain curve corresponding 
to the preferred configuration for a displacement-controlled problem. We see 
that a snap-through takes place at A = Ac, which leads to the localization 
(as represented by Solution 2). 




Fig. 10. The engineering stress-strain curve 
corresponds to the proffered configuration 

Once the localization happens, there is a high concentration of energy around 
the middle of the cylinder. It is known that if the strain energy density reaches 
a critical value there will be the material failure. The analytical results ob- 
tained here can be used to calculate the stored energy at any material point. 
The largest energy value is attained at {R, Z) = (a, ^) at which V = a (cf. 
Eqs. (61)-(63)). From Eqs. (40), (65), and (66), we can express the energy 
value at this point in terms of a, and the result is 

^-^ = -384FIT8D^(2A«'(-118 + 15a2 - 8192Dla^ 

-256(4 + 3a2)F2a + 16D3a(123 - 5a^ + 256(2 + a^)F2a)) 
+3{D2a*{-59 + lOa^ - AOdGDla^ - 512(1 + a^)F2a 
+4L)3a(246 - 150^ + 256(4 + 3a'^)F2a)) + 2(54K - 59a'^ 
+5a2a2 ^ ^^g^^ _ 5^2^^ ^ 4:096a'^{K - a{a - 2-f))Dl 

+256a{4:K + a(-(2 + a^)a + (4 + a'^)-f))F2 
+4L'3a(-2(118 + 5a^)K + a(246a - 5a'^a - 4927) 
+256a(2(-4 + a^)K + a((4 + a^)a - 8-f))F2)))) + 7a, 



The values of G^ corresponding to those values of A (in an increasing order) in 
Fig. 8 for preferred configurations are respectively Gm = 0.00722628, G^ = 
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(69) 



0.00742384, G„ = 0.00742989, and Gm = 0.00742989. 

Based on the maximum-distortion-energy theory (the Huber-Hencky-Von Mises 
theory; see Riley et al. 2007), there are two portions of the strain energy in- 
tensity: one for volume change and the other for shape change. In the present 
work, we consider an incompressible material, so there is no strain energy in- 
tensity corresponding to the volume change. Then the strain energy is only 
due to distortion. On the other hand, the strain energy intensity attains its 
maximum value at the material point (a, |). Thus, we can give the failure 
criterion 

Gm = Gf, (70) 



where Gf is the failure value of the strain energy intensity for a given material. 
Fracture will occur whenever the energy by Eq. (69) exceeds the limiting value 



8 Concluding Remarks and Future Tasks 



To study the localization phenomenon, a phenomenological approach is 
employed to formulate a three-dimensional boundary value problem with an 
incompressible hyperelastic constitutive law. A coupled series-asymptotic ex- 
pansion procedure is developed to solve the non-dimensionalized system of 
governing differential equations with given boundary data for a slender cylin- 
der subjected to axial tension. With the assumptions appropriate for the slen- 
der cylinder, analytical solutions have been obtained for the axisymmetric 
boundary value problem, which demonstrate the essential features of localiza- 
tion problems and are consistent with the experimental data available. Specif- 
ically, the width of localization zone depends on the material parameters, and 
it remains unchanged for the same material in the post-peak regime. Also, the 
snap-back and snap-through phenomena could be predicted with the analyti- 
cal approach, and a preferred configuration in the post-peak regime could be 
identified via an energy analysis. Due to the lack of three-dimensional ana- 
lytical results available in the open literature, the analytical work presented 
in this paper could complement the analytical, experimental and numerical 
efforts made by the research community for the localization problems over the 
last several decades. 

As indicated by Buehler et al. (2003), the hyperelasticity is crucial for un- 
derstanding and predicting the dynamics of brittle fracture. Especially, the 
effect of hyperelasticity is important for understanding the failure evolution 
in nanoscale materials. Since localization identifies the onset of material fail- 
ure, future work will focus on the identification of the parameters proposed in 
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the current phenomenological model, and on the hnkage between the contin- 
uum and fracture mechanics approaches, via muhiscale analysis. 
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